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ABSTRACT 


A modified form of the Lepping - Argentiero single spacecraft, shock 
normal determination procedure is presented. The modified method incorporates 
a simple predictor-corrector algorithm which allows a faster convergence rate 
and the use of average values of the parameters for the starting vector. 


INTRODUCTION 


A technique and associated computer program developed by Lepping and 
Argentiero (1971) to least-squares fit a sub-set of the Rankine-Hugoniot 
equations to shock-related plasma and magnetic field data, is modified in this 
note for the purpose of increasing the program's speed and efficiency and 
making unnecessary the need for multiple starting conditions in the iterative 
scheme. Originally the fitting scheme employed a standard Newton-Raphson 
numerical iterative procedure to solve simultaneously a set of eight nonlinear 
equations, denoted equations (36) in the original paper and equations ( 1 ) 
below. The iterative technique is described by Deutsch (19.65). When 
convergence did not result, one was required to repeat the procedure using, a 
different "starting vector" C Z Q ) , the eight components of which constituted 
the independent plasma and field variables to be adjusted. The simple 
modification suggested here, as well as making the use of more than one 
starting vector unnecessary, in general speeds up convergence by eliminating 
occasional large "overshoots" in the iterative procedure. 

Formulation 
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the input data for the eleven physical quantities of interest, M ( i ) data 
points each, usually arranged in the following order: 


'lx’ D 1y ’ D 1z’ 2x' 2y’ 2z’ 2y 
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(i.e., components of the magnetic field, plasma velocity differences and 
densities, where the subscripts 1 and 2 represent pre-shock and post-shock 
states, respectively), and. 

a i (i = 1, .. 11 ) 

are the "sigma noise parameters" used as weights in the least squares 
procedure and usually composed principally of the rms deviations of the 
measured physical quantities of interest (Y^'s), and finally 

X ± (i = 1, ... 8) 

are the independent variables (components of 2) to be solved for, physically 
representing the first eight quantities of Y^ . Note that X^ (aN^, X 1Q (=N 2 ) 
and X^ (= V 2x “ V lx^ are re -^ afced to (i = 1. ... 8) by the three lowest 
order Rankine-Hugoniot equations, which play the key role of providing a 
constraint on the nonlinear least-squares process [see equations (18), (19). 
and (20) in Lepping and Argentiero (1971)]. 


We denote Z (n = 1, n ) as an estimate of Z at iteration step n. 

n * max 

The components of Z q , the starting vector, are commonly the averages of the 
first eight components, Y^. Then at step n 




(3) 


When a "best estimate" of Z is obtained, the magnetic field components, i.e., 

A 

the first six components of that estimate of Z, are used in the magnetic 
coplanarity formula (Colburn and Sonett, 1966) to obtain a "best estimate" of 
the shock normal . 
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The original computer program is described in the Appendix of Lepping and 
Argentiero (1970). 

Modification of Scheme 


At each step n the quantity L, given by equation (2) and referred to as 
the least squares "loss function", can be evaluated. Our basic purpose is to 
minimize L(Z), as in equation (1). In order to avoid divergence of |AZ n l and 
L (2), we make a slight modification of the iteration procedure by allowing a 
choice of two branches in the fitting program, which are defined by the 
following : 


First branch - If at step n 


L„ < L . , then 
n n— 1 
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as in the original scheme. 


(4) 


Second branch - If at step n 


L > L , , then 
n n-1 ’ 


M = (L ,/L ) (Z - Z ,). 
n n-1 n n n-1 


(5) 


This continues step by step until n equals a pre-chosen integer, M say, or 
until 


( I A2 n |/ |f n - 1 |) < e (6) 

for some sufficiently small e > 0. Notice that the second branch differs from 
the unmodified first branch only in that the length of each step along 

unmodified) is shortened according to the ratio of the new to the old loss 
function. 

We have found in numerous cases that this simple modification has 
decreased, and sometimes significantly, the time needed to run the fitting 
program and/or the operator's efforts in finding a proper starting vector. In 
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extreme cases, for a given 2 , the original program would occasionally wander 
around the neighborhood of L^ n with there being no apparent hope of 
satisfying (6) for a reasonable z, especially when the gradient of L(2) was 
shallow near L m ^ n , i.e., for poorly conditioned cases. This modification 
should eliminate that problem. Although not a unique solution to the 
convergence problem, this type of modification is attractively simple. 

Example 


Here we show the benefit of using the modified program by comparing its 
results with those of the unmodified program for a somewhat poorly behaved 
case; we have encountered much worse cases but they are not typical and 
therefore not good examples. 

We examine a shock observed on Pioneer 6 at 2058 U.T. on March 22, 1966, 
which was first studied by Chao (1970); he graphically displays and discusses 
the data, which are from the Goddard magnetometer and the MIT plasma probe 
onboard the spacecraft. The right side ("DATA") of Table 1 shows the 
tabulated data points constituting 12.5 minutes of upstream field data, 14 
minutes of downstream field data (partially decimated), and ±22.5 minutes of 
all plasma quantities. The analysis intervals were chosen on the basis of 
apparent steadiness of the data. All data are given in an R-T-N coordinate 
system, centered at the spacecraft, where ft is positive radially away from the 
sun, T is perpendicular to ft , parallel to the ecliptic plane and positive in 
the direction of the earth's motion about the sun, and ft = ft x T. On the left 
side of the table are best estimate (B.E.) values from both the modified and 
unmodified programs, straight averages (AVG) of the quantities, and the 
associated sigma noise parameters (SIG), which in this case were simply rms 
deviations of the data values. The modified program converged to the B.E. 
value in 3 iterations where the AVG's were used for the components of the 
starting vector. By contrast the unmodified program required 11 steps when 
again the AVG starting vector was used. In the unmodified program the loss 
function L had the undersirable feature of oscillating in value, as shown in 
Table 2; the table also shows the monotonically decreasing values of L that 
the modified program produced, as expected. The table also gives a related 
quality factor defined below the table. We stress the case where AVG's were 



used for the starting vector, Zo, since that is obviously the easiest and most 
commonly used choice, as stated above. However, for this shock several other 
apparently reasonable Z Q ’s were attempted using the unmodified version of the 
program, and seme did not give convergence after 15 interation steps; the 
program was set to stop at the 15th step. 

It must also be emphasized that the selection of an optimum fitting 
coordinate reference frame is important in assuring proper and speedy 
convergence. The implementation of the modified and unmodified procedures 
have included interactive facilities to conveniently rotate the data into one 
of three orthogonal alignments such as to assure the optimum selection of 
dependent and independent variables in the fitting procedure. In general, the 
quantities exhibiting the greatest variability or uncertainty should be 
selected as the dependent variables to be least-squares estimated. 

The best estimate shock normal, based on the B.E. values in Table 1, 
is 


fl B.E. S (n R' n T’ V S (0 - 94 ’ -°* 15 » 0-3D. 

and the associated 95* certainty error cone half-angle was 7.7° (Lepping and 
Argentiero, 1971). [If AVG value fields had been used, the error cone 
half-angle would have been 3.3 times larger!]. By comparison Chao’s estimate 
of the shock normal was 

A 

n CHA0 = (0 * 84 ’ ” 0 * 10 ’ 0,54) ’ 

differing by 13.6° from our estimated direction. The shock was a rather 
typical oblique one (at 1 AU), whose B.E. normal was 59° from the upstream 
field direction B^^gg. and whose magnitude ratio across the shock was 

!b( 2 ) be 1 /|b( 1 ) se ! = 1,8 °* 

We stress how different the velocity difference vector W = (WR, WT, WN) 
is between the B.E. and AVG values, especially in the WN component, as Table 1 
shows. The angle between Wgg. and W AV£ is 17.1°. If one were to use "velocity 
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coplanarity" ( Abrahara-Shrauner , 1972) and AVG values, the resulting shock 
normal would be 

"vel S S AVG /|S fLVG l * (0 - 69 ' °- 70) ’ 

which differs by 26.4° from our best estimate normal; the author warns that 
this method is an approximation. Our added warning is that using average 
values is often inadequate. For example, we see that by using the B.E. W and 
velocity coplanarity we obtain 

n vel 5 i? BE /| ^3E l = (0 * 86 ’ “°* 23, 0<46) ’ 

A 

which is only 9.8° from n 3 £. Obviously in all of the above we have made the 
tacit assumption, argued in Lepping and Argentiero (1971), that is indeed 
generally the best estimated shock normal from data from a single spacecraft. 

Many such examples can be found to show the desirability of using the 
modified version of the ’best-fit’ technique. 
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TABLE 2 


Loss Function Values 


Iteration Step L-unmodified (Q*) 


L-modified (Q*) 


0 

{Starting vector 

1,095 

(.295) 

455 

(.457) 

1 


19,360 

(.070) 

174 

(.739) 

2 


2,985 

(.178) 

111 

( . 925 ) 

3 


2,749 

(.186) 

108 

(.938) 

4 


7,336 

(.114) 



5 


805 

(.344) 



6 


5,013 

(. 138) 



7 


1,156 

(.287) 



8 


1,639 

(.241 ) 



9 


127 

( . 865 ) 



10 


139 

( . 826 ) 



11 


108 

(.938) 



*Q s 

Quality = (N t /L) 1/ ? 

where is 

the total 

number of data points. 

95 in 

this 

case. Either two sucessive Q's 

of 0.85 or 

larger , or Q of 0. 90 

or 


larger, is usually considered a successful convergence, provided reasonable 
"sigma parameters" were used. 
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